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ABSTRACT 

We perform N-body simulations of the dynamical evolution of a dwarf galaxy 
falling into the Milky Way galaxy in order to understand the formation scenario of 
the peculiar globular cluster uj Centauri. We use self-consistent models of the bulge and 
the disc of the Milky Way, as well as of the dwarf galaxy, and explore a range of dwarf 
models with different density distributions. Namely, we use King (1966) and Hernquist 
(1990) density profiles to model the density distribution in the dwarf. The central 
region of our King model has a density profile approximately cx while that of the 
Hernquist model is cx r^^. The difference in the dwarf's density distributions leads 
to distinct evolutionary scenarios. The King model dwarf loses its mass exponentially 
as a function of apocentric distance, with the mass-loss rate depending on the initial 
mass and size of the dwarf. Regardless of the initial mass and size, the King model 
dwarf remains more massive than 10^ Mq after a few Gyr of evolution. The Hernquist 
model dwarf experiences an accelerated mass loss, and the mass of the remnant falls 
below 10^ Mq within a few Gyr. By exploring an appropriate set of parameters, we 
find a Hernquist model that can attain the mass and orbital characteristics of uj Cen 
after a few Gyr. 

Key words: galaxies: dynamics and kinematics - galaxies: interactions - globular 
clusters: individual: lu Centauri - methods: A'^-body simulations 



1 INTRODUCTION 

The accretion and subsequent tidal disruption of dwarf 
galaxies are important agents in the formation and evolu- 
tion of the Milky Way (Searle & Zinn 1978). Such accretion 
events would presumably lead to formation of tidal streams 
in the halo of the Milky Way associated with the disrupted 
dwarf galaxies. The discovery of the tidal stream related to 
the Sagittarius dwarf galaxy (Ibata et al. 2001) indeed sup- 
ports the idea that halos of galaxies can have a number of 
streams maintaining their coherence during a few Gyrs. 

There is another class of halo objects which might be 
associated with past accretion events. During the accretion 
of nucleated dwarf elliptical galaxies, the compact stellar nu- 
clei may survive the disruptive event and continue orbiting 
until the present time. Freeman & Bland-Hawthorn (2002) 
even speculate that the globular clusters in the Milky Way 
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could be the stripped relics of ancient protogalactic stellar 
systems. To identify such objects would be of great impor- 
tance in understanding the nature of such building blocks 
of galaxies. In this paper we focus on one object which may 
be the relic of a past accretion event - the globular cluster 
LU Centauri. uj Centauri is one of the most peculiar globular 
clusters in the Milky Way: with a mass of 5 x 10® Mq , it is 
the most massive globular cluster. It has an unusual flatness 
with an e llipticity of abo ut 0.12 associated with fast internal 
rotation l)Freemanll200l . It orbits in a retrograde direction 
relative to the disc with an apocentric radius of about 6 kpc, 
a pericentric radius of 1.2 kpc, and a ma ximum height above 
the plane of 1 kpc (iDinescu et al.|[l999ll . 

What really distinguishes uj Cen from other Galactic 
globular clusters, is its chemical abundance pattern. It is the 
only known globular cluster that has well-defined signatures 
of self-enrichment. Thus, according to Smith et al. 2000 and 
Vanture et al. 2002, the abundance pattern observed within 
LO Cen can be understood as a combination of Type II super- 
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novae and, more importantly, AGB-star enrichment during 
a period of about 3 Gyr. Gnedin et al. (2002) have demon- 
strated that uj Cen on its present orbit and with its present 
mass can not retain heavy elements dispersed in the AGB 
phase of stellar evolution. The passages through the disk 
would sweep out the intracluster gas and the cluster, once 
established on this orbit, would not be able to produce stars 
with enhanced s-process elements. Thus, the cluster should 
have evolved chemically before settling on its present orbit, 
or, in other words, it was accreted by the Galaxy. Another 
possibility to explain the chemical peculiarities of u Cen is 
the merger of two or more chemically distinct systems. Van- 
ture et al. (2002) conclude however, that in this scenario, it 
would be difficult to explain an increase in the of ratio of 
heavy s-process elements to iron with metallicity. 

A plausible explanation is that uj Cen is the nucleus 
of a dwarf galaxy captured and disrupted by the the Milky 
Way (e.g.. Freeman 1993). In this scenario, a dwarf sinks to 
the center of the Galaxy due to dynamical friction, simul- 
taneously losing mass in the Galactic tidal field. All of the 
dwarf's mass is eventually stripped but the nucleus, with a 
mass of 5 X 10® Mq , is left on the current orbit of uj Cen. 
The current angular momentum of uj Cen is thus inherited 
from the dwarf galaxy, on a retrograde orbit. The complex 
chemical composition and extended star formation history 
of u! Cen can be explained in this picture by the chemical 
evolution of the nucleus in a deep potential well of a dwarf, 
such that the nucleus would be able to retain stellar ejecta 
during subsequent activity of supernovae and AGB stars. 

Even though a capture scenario is seen as a plausible 
explanation for the peculiar properties of uj Cen, it is unclear 
whether it is possible to realize an orbit and mass evolution 
of the dwarf that can reproduce the current parameters of cu 
Cen. Zhao (2002) studied the orbital decay of a dwarf that 
may lead to the current position of cu Cen by launching a 
dwarf on an orbit that started 50 kpc away from the Galactic 
center. Using a semi-analytical model, he concluded that a 
progenitor of uj Cen can not decay to its present orbit. He 
found that strong tidal shocks quickly reduce the mass of 
the dwarf so that dynamical friction becomes too weak to 
drag the remnant to the inner regions of the Galaxy. Zhao 
(2002) concluded that the only possibility to explain uj Cen 
- phenomenon is that its progenitor was born only 15 kpc 
away from the Galactic center. 

Zhao's analysis, however, was simplified in a number of 
ways. He assumed that the mass outside the tidal radius 
is instantaneously stripped, and the contribution from the 
bulge and disc components to the dynamical drag was not 
taken into account. In this study, we perform a series of 
numerical simulations, using self-consistent bulge and disc 
models for the Milky Way, and self-consistent dwarf mod- 
els that have different density profiles. As a result, we find a 
dwarf-capture scenario that successfully produces an uj Cen- 
like object from a normal dwarf galaxy. Th is result is also 
repor ted in more detail in a separate paper jTsuchiva et alJ 
l2003r . In this paper, we explore the model parameters and 
focus on the physical processes of the dynamical evolution of 
such infalling dwarfs making this study beneficial in under- 
standing the properties of merging parent-satellite galaxy 
systems. 

The numerical realization of this scenario is a compli- 
cated problem because of the large range of masses and sizes 



of the interacting components, and due to the drastic change 
of the mass of the sinking dwarf. We explain our numerical 
methods in Sect. |21 where we describe the Galaxy and the 
dwarf models. The results of the simulations are presented 
in Sec. 13 Some analytical considerations for our results are 
given in Sect. 4. Finally we formulate our conclusions in 
Sect.ini 



2 MODELS AND NUMERICAL METHODS 

2.1 Af-body simulations 

In our Ai'-body simulations, we neglect gaseous effects such 
as pressure and viscosity. We do not include star formation, 
and do not distinguish between visible and dark matter in 
the models. Gravitating matter is represented by the soft- 
ened particles, and their mutual gravity is calculated by a 
hierarchical tree algorithm (see, for example Dubinski ( 199^) 
for a review). The tolerance parameter which is used to or- 
ganize the hierarchical boxes, is 0.7. The softening length is 
0.01 of the disc scale length, which corresponds to 35 pc in 
our Milky Way model. We use the same value of the soften- 
ing length for all particles irrespective of their mass. 

The integration time step is 1/32 of the dynamical time 
which is determined as the quarter of the orbital period of 
a particle at the center of the Milky Way. This gives a value 
of 4.38 X 10^ years. The central density of our dwarf is an 
order of magnitude larger than that of the Milky Way model. 
However, the orbiting time of the stars in the dwarf is longer 
compared to the time step. 

An oscillation time of the dwarf in the potential of the 
Milky Way is about 0.5 Gyr at the beginning of an infalling 
process, and about 0.1 Gyr at later stages of evolution. Since 
the plane of the orbit of the dwarf has a low inclination with 
respect to the disc, a characteristic time of the gravitational 
shocking from the disc is comparable to the dwarf's orbit- 
ing time. The bulge shocking time determined as a time of 
crossing of the bulge by a moving satellite is much shorter 
(about 3 X 10® years) compared to the gravitational shoking 
time from the disk, making bulge shocking much stronger. 

2.2 The Milky Way Model 

An accurate model of the Milky Way is important to 
make a quantitative estimate of the dwarf evolution. We 
construct the Milky Way equilibrium model following 
Kuiikcn and Dubinski ( 1995). The model provides wide flex- 
ibility in fltting observational constraints. Briefly, it consists 
of a nearly spherical bulge, an exponential disc, and a halo 
which resembles a lowered isothermal sphere at larger scales. 
Each component is expressed by a distribution function, so 
that the positions and the velocities of A'^-body particles can 
be directly sampled from the distribution functions. Since 
the model is close to equilibrium, there is no need to "relax" 
it before an actual simulation. 

The parameters of the model are chosen to reproduce 
the properties of the real Milky Way galaxy. We choose the 
solar radius, the disc exponential scale length (Rd), and the 
disc vertical scale height (za) to be 8 kpc, 3.5 kpc, and 245pc 
respectively. The circular velocity of the disc at the solar 
radius is taken to be 220 km/s. The total surface density 
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Figure 1. The rotation curve of our Milky Way model. The con- 
tributions from the bulge, disc and halo to the circular velocity 
are shown by dotted, short-dashed, and long-dashed curves, re- 
spectively. 
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Figure 2. Artificial heating of the Milky Way disc model due to 
two-body relaxation. The number of particles in the bulge and 
disc components is 10,000 and 70,000 respectively, and the halo 
is treated as a fixed potential. The softening length in the tree 
code is 35 pc. 



within 1.1 kpc of the disc plane, and the contribution of the 
disc material to it are 69.8 Mq pc~^, and 45.5 Mq pc~'^ re- 
spectively, which agrees with the obser vational constraints of 
71 ±6 Mq pc"^ and 48 ±9 Mq pc"^ jKuiiken and Gilmord 
I1991I) . With a cut-off radius of 28 kpc, the disc total mass is 
5 X 10^" Mq . Toomre's Q value of the disc is about 1.9 at a 
distance of 2.5 Rd from the center. The bulge mass and the 
cut-off radius of the bulge are taken to be 0.75 x 10^" Mq 
and 2.38 kpc respectively. 

The halo density distribution can be well ap- 
proxima ted by a spherical lowered isothermal profile 
tBinncv and Tremaine 1987.') with a central potential 
^'(0)/(T^ = 8 except for the central regions which are de- 
formed owing to the presence of the bulge and the disc. The 
halo mass inside 50 kpc and 170 kpc radii is 4.9 x 10^^ Mq 
and 8.6 x 10^^ Mq respectively. This is within the obser- 
vationally inferred values 5.4l3;g x lO" Mq and x 
10^^ Mq found by Wilkinson and Evans (1999). The contri- 
bution of each component to the circular velocity is shown 
in Fig. □ 

This model is similar to model S in iTsuchival ll2002l) . 

An interested reader can find in that paper a more detailed 
description of the model and its stability analysis. The bulge 
and the disc in our sumulations are represented by 10,000 
and 70,000 equal mass particles. With this rather low num- 
ber, the disc thickness grows due to two-body relaxation. 
Fig. |5]shows the change in the disc thickness in equilibrium. 
During 6 Gyrs, the disc thickness increases up to 500 pc near 
_R = 10 kpc. This is about twice as large as the present day 
thin disc, but smaller than the scale height of the thick disc. 
The bump in the disc thickness around R — 3 kpc is caused 
by a weak bar developing in the disc. A further discussion 
of this effect can be found in Tsuchiya (2002). 

In this study, the halo is treated as a fixed potential 
since we are mainly interested in the effects of the disc and 
bulge shocking and dragging of the dwarf galaxy. We ne- 
glect dynamical friction from the halo in our models, which 
increases the orbital sinking time (see also Section 5.1). If 
we treat the halo as an ensemble of A'^ particles, additional 



heating causes disc thickening. To prevent considerable disc 
heating, the number of particles in the halo should be sev- 
eral times larger than that in the disc. This implies severe 
limitations on our computational capabilities. Nevertheless 
we run one simulation with a "living" halo to estimate the 
contribution of the halo to the dwarf's dynamics. The halo 
dynamical friction decreases the dwarf's orbital decay time 
and increases its mass loss, but does not change the evolu- 
tion history in the mass-apocenter radius plot. We discuss 
this issue in section 5.1. 



2.3 Dwarf Models 

We assume that a dwarf galaxy captured by the Milky Way 
was the hypothetical progenitor of uj Cen. Presumably, a 
high-density nucleus of the dwarf survived the disruptive 
capture event and ended as u Cen. We wish thus to ex- 
plore dwarf models with centrally concentrated density pro- 
files. We exami ne two models: the lowered isothermal or 
the King mo del l|B irmgY mkI Tre mainell98^ . and the Hern- 
quist model llHerna uist"l990^. Both models are spherically 
symmetric. The kinematics of the compact dwarf spheroidal 
galaxies does not reveal a significant rotation of these sys- 
tems (Binney & Merrifield 1998). We assume therefore that 
both of our models have no rotation. The King model is 
defined by the distribution function 

Ms) J ^>0; 

[ e < 0, 

where e = + is a specific relative energy of the 

particles defined so that the energy is zero at the surface 
with zero volume density. The value of the dwarf potential 
at the center, '5(0), determines how centrally concentrated 
the distribution is. We adopt the value *(0)/(t(0)^ = 12. 

The model has two charact eristic scale lengths: the 
King radius rx = •\/9(Tg/47rGpo, and the tidal radius n. 
For ^'(0)/(To — 12, the ratio between the tidal and the 
King radii, is rt/vK = 548. These two scale lengths are 
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Figure 3. Density distribution (upper panel) and velocity dis- 
persion (lower panel) in tlie King model witli 'I'(0)/(Tq = 12. The 
total mass and the half-mass radius of the model are 8 X 10^ Mq 
and 4 kpc respectively (model K5). The distribution is taken from 
our Af-body realization with N = 50, 000 particles averaged over 
spherical shells containing 100 particles for the density distribu- 
tion, and 1000 particles for the rms velocity distribution. The 
dashed line in the upper panel shows the density distribution in 
the model after 5.5 Gyr of evolution in isolation. For comparison, 
we also show the density distribution (dotted line). In the 
lower panel, the dashed and dotted lines show the radial and the 
tangential velocity dispersion profiles respectively, after 5.5 Gyr 
evolution of the model in isolation. 



not convenient however, for the comparison with the Hern- 
quist models, so we label the distributions by their half- 
mass radii related to the King radius and the tidal radius 
as ri/2 = 0.164rt = 89.9rK- Figure |3 shows the density and 
velocity profiles in our King model with total mass 8 x lO' 
Mq and half-mass radius ~ 4 kpc sampled by 50,000 
particles. Since the model has a deep central potential, a few 
hundred particles compose the central homogeneous core of 
the dwarf. 

To study the numerical stability of the model, we have 
followed its evolution in isolation for about 5.5 Gyr. The 
final distributions are shown in Fig. |3 The central density 
decreases by a factor of 2. The change in density profile is 
noticeable within r < 0.1 kpc of the dwarf center which con- 
tains mass of about 10* Mq . The central velocity disper- 
sion decreases by about 25%, while there is no considerable 



change in the velocity dispersion profile in the outer regions 
of the dwarf, i. e., at r > 0.5 kpc. 

The Hernquist model (Hernquist 1990) has a density 
profile given by the equation 



p{r) = 



2irr{r + ro) 



(2) 



In central regions, the density profile has a cusp of r~^, 
while in the outer regions the density decreases as r~*. The 
radius ro separating the two regions is related to the half 
mass radius as: 

ri/2 = (1 -f V2)ro. (3) 

The slope of the density gradient in the central region 
is thus shallower than that in the King models. 

A more noticeable difi'erence between the two models is 
in their velocity dispersion profiles. While the King models 
have constant velocity dispersion in the center, in the Hern- 
quist models the velocity dispersion decre ases near the cen- 
ter as (T^ oc r/ro ln(ro/r) iHernauislll990t) . Figure^Jshows a 
realization of the Hernquist model which has a total mass of 
Mtot = 8 X 10^ Mq and half- mass radius of ri/2 = 2\/2 kpc. 
This Figure also shows the density and the velocity disper- 
sion profiles after 5.5 Gyr evolution in isolation. The density 
in the central region, within the softening length, decreases 
in the same way as for the King model. However, the velocity 
dispersion profiles remain close to the initial disrtributions. 

Evolution of the dwarf's central regions is caused by the 
softening of the gravity forces. As can be seen from Figure 
3, the velocity dispersion profile in the King models fiattens 
toward the center, while the volume density grows faster 
than r~^. To keep the central regions in equilibrium, the 
self-gravity force in the King models should grow towards 
the center. The velocity dispersion in equilibruim Hernquist 
models decreases toward the center( Figure 4), which results 
in a constant self-gravity force in the central regions. The 
softening length of 35 pc is comparable to the King radius 
of 44.5 pc. Thus, the gravity softening decreases the gravity 
force in the core regions, and causes the core to expand re- 
ducing its density and velocity dispersion. Central regions of 
our numerical King equilibrium models with softened grav- 
ity are more departed from equilibrium than those of Hern- 
quist models, which results in more noticeable evolution of 
central regions of the King models. 

Decreasing the density at the center of the dwarf can 
potentially cause its more rapid disruption. We find how- 
ever, that the mass of the King models remains larger than 
10* Mq during the simulations and a structural change at 
the dwarf center does not affect mass stripping in the simual- 
tions. Our test simulations thus show that isolated dwarf 
models are fairly stable at scales larger than a few softening 
lengths. In other words, artificial effects such as the soft- 
ening and the two-body relaxation do not affect the mass 
stripping history unless the mass of the remnant falls below 



10^ Ml 







2.4 Multi-Mass A/^-Body Model of the Dwarf 

With our choice of the softening length, the numerical reso- 
lution near the dwarf's center is about 10^ Mq . In our King 
models, we use 50,000 particles to represent a 8 x 10^ Mq 
dwarf, and a lO'^ Mq object would consist of a few tens of 
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Figure 4. Density (upper panel) and velocity dispersion (lower 
panel) distributions for the Hernquist model with total mass of 
8 X 10^ Mq and the half mass radius of 2.83 kpc (model H5) 
taken from multi-mass Af-body realization with A'^ = 140, 000 par- 
ticles. The density distribution is determined from averages over 
spherical shells that contain 100 particles each, while the veloc- 
ity dispersion distribution from shells that contain 1000 particles. 
The dashed line in the density plot represents the distribution 
of the model after 5.5 Gyr of evolution in isolation. The lines 
in the velocity dispersion plot show the radial (dased) and tan- 
gential (dotted) velocity dispersions after 5.5 Gyr of evolution in 
isolation. 



particles. To improve mass resolution, we increase the num- 
ber of particles close to the dwarf center. We find that the 
Hernquist models undergo strong mass loss, and the mass 
of the remnant becomes less than lO'^ Mq during a few Gyr 
while the King model dwarfs stop losing their mass above 
10* Mq . To avoid effects caused by a small number of par- 
ticles at the late stages of evolution, a remnant should be 
composed of more than a thousand particles. To increase the 
mass resolution of the Herquist model near the center, we 
use a mulyi-mass N-body model of the dwarf. We divide the 
particles in the model into three groups differing by their 
energy: 



fH{E)dxdv = 0.01, 



E<Ei 



fH{E)dxdv = 0.9, 



E2<E 



fH{E)dxdv = 0.09, 



(4) 



I Ei<E<E2 



where E — + $(a;) is a specific energy of the in- 
dividual particles. The lowest energy group is sampled by 
the particles with mass mi = 2 x lO"'^ A/tot, and with the 
number of the particles in the group equal to 50,000. The 
middle and the high energy groups are sampled by the parti- 
cles with masses 7712 = 2 x 10~^Aftot and ms = 2 x 10~^Mtot 
respectively, and with the number of particles in each group 
of 45,000. Since the light particles are progressively con- 
centrated near the center, the mass resolution of the dwarf 
central regions becomes about a hundred times better. The 
central lO'^ Mq for example is composed now of more than 
50,000 particles. 

A potential problem of N-body models with different 
masses of particles is an artificial heating of light particles. 
In our model, light particles are concentrated toward the 
dwarf's center, while heavy particles are distributed at the 
periphery of the dwarf. On the contrary, in a dwarf close to 
equilibrium, massive particles sink toward its center, while 
light particles are scattered away. In other words, there is a 
net energy fiow from the massive particles to the light ones 
in our model. We have checked whether this process is sig- 
nificant in our simulations. Fig. |3 shows energy exchange 
between particles of different masses as a function of time 
for the equilibrium multi-mass Hernquist model evolving in 
isolation. Solid, dashed and dotted lines in this Figure stand 
for the low, medium and high mass particles. The left panel 
in Fig. shows the total energy of each group of particles 
as a function of time, while the right panel shows the spe- 
cific energy, averaged over the particles in each component. 
Ninety percent of the mass of the dwarf in our model is 
built with massive particles. Therefore the specific energy of 
heavy particles does not change significantly, while the total 
energy does. 

There is no noticeable relaxation between the low and 
intermediate mass particles. The global energy flow is domi- 
nated by the energy transfer from heavy to the intermediate 
and low-mass particles with the dominant energy exchange 
being between the massive and intermediate mass ones. The 
light particles also gain energy as a result of relaxation. This 
energy gain, although, does not cause a significant change 
in the equilibrium density ditribution as can be seen in Fig- 
ure 4. Furthermore, heating of the light particles is even less 
important in our simulations since most of the heavy parti- 
cles are stripped away from the dwarf at early stages of the 
infalling process. 



2.5 Bound mass 

The bound mass of a dwarf moving in an inhomogeneous 
background field is a rather ambiguous quantity. One way 
to determine the dwarf's bound mass is to calculate the 
total mass of the dwarf particles that have negative energy, 
ignoring the external gravitational field of the Milky Way: 

-^SG = ^'rri'ivf + ^dwarf (a;0 < 0. (5) 

Another way is to calculate the mass inside the tidal 
radius where the tidal force from the Milky Way is com- 
parable to the self-gravity of the dwarf. This would be an 
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Figure 5. Energy change in the three different mass components 
of the multi-mass Hernquist dwarf. Panel (a) shows the total 
energy of each component, while panel (b) shows the specific 
energy averaged over the particles in each component. The solid, 
dashed, and doted lines stand for the light, medium, and heavy 
mass components respectively. 



appropriate bound mass definition if the dwarf is subject 
to a strong tidal field. To simplify the tidal radius calcu- 
lations, we assume that the Milky Way is spherically sym- 
metric with a homogeneous core of 1-kpc radius, and has 
a flat rotation curve of 220 km/s outside the core region. 
We use the expression for the tidal radius from Binney and 
Tremaine (1987) ignoring the deviation of the Milky Way 
from spherical symmetry: 



f Mdwarf(n) \ 

3MMw(d) ) 



1/3 



(6) 



Here, d is the distance between the center of the dwarf (de- 
fined at the location of the minimum of the dwarf's po- 
tential) and the center of the Milky Way, Mdwarf ('"Oi s^nd 
Myrw{d) are the mass of the dwarf within the tidal radius, 
and the mass of the Milky Way within radius d. Under these 
assumptions, the dwarf tidal mass is determined by the ex- 
pression: 



Md 



tin) 



3.27 X 10' 
3.27 X 10^ 



Vrot 



\ 220kms 

/ Vrot 



Mo 



.kpc^ 



d< 1 kpc 
d > 1 kpc 



(d/1 kpc)2 V220kms-i; Vkpc' 

As a practical compromise, we also calculate gravita- 
tionally bound mass of the dwarf within its tidal radius. 
Fig. Qshows the bound mass of an infalling dwarf as a func- 
tion of time defined in this way, as well as the bound mass 
calculated with the other two definitions. The bound mass 
determined within tidal radius of the dwarf does not have 
peculiar features observed in the mass profile calculated for 
the mass of the remnant bound by its self-gravity. 



3 RESULTS 

This section presents the results of the numerical simula- 
tions of the infall of a dwarf that has King and Hernquist 
density profiles. To match approximately the current orbit 
of uj Cen, the initial position and velocity of the dwarfs 
have been choosen to be {X, Y, Z) — (50, 0, 30) kpc, and 
{Vx, Vy, Vz) = (0, -20, 0) km/s. Here X, Y, Z are the Carte- 
sian coordinates associated with the Galaxy, X is positive 
away from the Galactic center, Y is positive in the direction 
of the Galactic rotation, and Z is positive toward the North 
Galactic pole. The total masses and the half mass radii of the 
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Figure 6. Time evolution of the orbital radius of the dwarf model 
K5. 



examined models are listed in Table 1 for the King models, 
and in Table 2 for the Hernquist models. These parameters 
control the overall evolution of the dwarf. Massive dwarfs 
sink faster in the potential well of the Milky Way, while 
dwarfs with the large half-mass radii, i.e. a lower volume 
density experience more rapid mass loss. 

The King models (Kl, K4), (K2, K5) and (K3, K6) 
ifsVed in Table 1 are related to each other as 



A/d 



f/r-i 



/2 



: const. 



(8) 



With the density profile of King models close to p oc r~ , the 
pairs of models have the same radial density distributions, 
and implicitly same central densities. Thus, pair (Kl, K4) 
has the highest central density, while (K3, K6) has the lowest 
central density. The Hernquist models (HI, H4, H7), (H2, 
H5, H8) and (H3, H6, H9) listed in Table 2, are related as 



Mdwarf/''l/2 = const. 



(9) 



Similarly to the King models, these models have the same 
density distribution and therefore central density within 
each group. Group (HI, H4, H7) has the highest central 
density, while group (H3, H6, H9) has the lowest. 

3.1 King Models 

All the King models listed in Table 1 are composed of 
50,000 identical particles and have central potential depth 
of *(0)/crg = 12. 

Figure 6 shows the time evolution of the dwarf orbit for 
model K5 taken as an example. This plot shows the distance 
of the dwarf from the center of the Milky Way as a function 
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Figure 7. The mass of the dwarf model K5 as a function of time. 
The soUd and dashed Hnes show the mass of the dwarf calculated 
as gravitationally bound and tidal mass respectively. The dotted 
line is the bound mass determined as gravitationally bound mass 
of the dwarf within its tidal radius. 



of time with the position of the dwarf determined as the 
point where the dwarf potential reaches its minimum. 

The dwarf has a highly eccentric orbit which decays due 
to the interaction with the disk and the bulge of the Milky 
Way. The distance minima in Figure 6 do not correspond 
to the exact pericentric passages of the dwarf because of 
the rather low time resolution in the plot, that is about 700 
Myr. We estimate the typical dwarf's pericentric distances 
of about 1 kpc. 

Fig. |7| shows the mass loss history for the dwarf model 
K5. We use three different definitions of the dwarf's bound 
mass. The dashed line shows the mass of the dwarf within its 
tidal radius. The dotted line shows the behavior of the grav- 
itationally bound mass of the dwarf within its tidal radius. 
For comparison, we plot also the mass of the gravitationally 
bound particles of the whole dwarf (solid line). The mass 
of the remnant determined by the gravitationally bound 
particles is larger than the mass of the remnant within its 
tidal radius. Until t ~ 1.4 Gyr, the mass of gravitationally 
bound particles of the dwarf decreases monotonically, chang- 
ing later on into an oscillatory regime. By defining the rem- 
nant mass as the mass of the dwarf's gravitationally bound 
particles, we neglect the gravity of the Galaxy. With this 
definition, the particles with sufiiciently small velocities will 
be gravitationaly bound even if they are located at large dis- 
tances from the dwarf center. Close to apocenter, the dwarf 
and the stripped particles slow down, which results in an 
increase of the dwarf bound mass. We find this definition 
inappropriate for our problem. 

The mass of the remnant determined by its tidal radius 
changes mainly during its pericentric passages and remains 
nearly constant for the rest of the orbiting time. A simple 
approximation for the tidal radius estimate is invalid how- 
ever close to pericenter, where the dwarf's distance from the 
Galactic center is comparable to the size of the dwarf itself. 
The bound mass determined as a mass of gravitationally 
bound particles within the tidal radius ofthe dwarf is thus 
the most plausible definition of the gravitationally bound 
mass for our problem. 

The orbital sinking of the dwarf is caused by dynamical 



10 



10 



10=^ 



« 10*^ 



10' 



K3 ^ 

K6 - - « 


) — 


! 








" 


....... -I 


K4* 















Iomega Gen 







10° 

10 20 30 40 50 60 
Radius (l^pc) 

Figure 8. Evolution of the mass as a function of apocentric dis- 
tance for the King models. The labels in Figure match the models 
hsted in Table 



friction, or a back reaction from the disc and bulge den- 
sity perturbations determined by the dwarf's gravitaty. In 
an ideal case, the drag force is proportional to the square of 
the dwarf mass (Chandrasekhar 1943). As the dwarf loses 
its mass, the orbital decay rate decreases. The history of the 
orbital decay depends thus crucially on the mass loss rate of 
the dwarf, which in turn depends on dwarf's internal den- 
sity distribution. We study these complex connections in an 
evolutionary history of the dwarf by a systematic parameter 
survey of our models. 

An illuminating diagram is Figure 8 which shows the 
bound mass as a function of apocentric distance of the dwarf. 
The evolution history of King dwarf models can be classi- 
fied into the three groups depending on their initial central 
density. Models in each pair, e.g, (Kl, K4), (K2, K5), and 
(K3, K6) reach nearly the same bound mass and apocentric 
distance in their evolutionary courses. These model pairs 
have the same density distributions, and once the outer 
parts of the dwarf are stripped, the models in each pairs 
have nearly identical remnants. The bound mass of the King 
model dwarfs decreases exponentially with apocentric ra- 
dius. An extrapolation of the evolution of model K2 could 
lead in principle to a remnant with the mass and the orbital 
parameters close to those of uj Cen. Ifowever, the evolution 
time to reach to Gen's parameters is larger than the age of 
the Galaxy, and uj Gen's origin can not be explained by this 
model. 

Figure 8 also illustrates a common evolution of the mod- 
els that have the same central densities. Models K2 and K5, 
for example, have close masses by the second pericentric pas- 
sages. Models K3 and K6 have close bound masses even at 
the beginning of the simulations since this pair of models 
has the lowest central density, and the tidal radius for these 
models is already small at J? ~ 60 kpc. The additional initial 
mass in model K6 compared to the model K3 is distributed 
outside the tidal radius (as defined in Section 2.5), and is 
stripped during the very first passage by the Galactic center. 

Figure 9 illustrates the mass evolution of the King dwarf 
models. By 6 Gyr, the most the dwarf models reach mass of 
about few X 10* Mq, with the mass loss comparable with that 
found in Zhao' model (Zhao 2002). Although the models K3 
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Figure 10. Evolution of the angular momentum as a function of 
the decaying apocentric radius of the King dwarfs. 



and K6 have higher mass loss rates and reach 10* Mq during 
4 Gyr of evolution. 

A similar evolutionary pattern among different King 
models can be seen also in dependence of the specific an- 
gular momentum L^. Fig. IIUI shows Lz with the apocentric 
radius. The positive Lz is for retrograde orbits to the disc 
rotation. At the beginning, we set the same orbital param- 
eters for all models so that all models have the same initial 



3.2 Hernquist Models 

To examine the dynamics of the Hernquist model dwarfs, 
we build nine models listed in Table 2. As it was mentioned 
before, there are three groups of models namely (HI, H4, 
H7), (H2, H5, H8), and (H3, H6, H9) which differ in their 
half-mass radii and mass but have the same radial density 
distributions (i. e., central densities). 

Figure 11 shows mass evolutionary tracks for the nine 
Hernquist models. In this Figure, the bound mass of the 
remnant is plotted as a function of its apocentric distance. 
Similarly to the King models, a common behavior of the 
Herquist models which have the central density is seen in 
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Figure 11. Once the tidal interaction strips the outer lay- 
ers of a dwarf belonging to the same group, the remnants 
become nearly identical and have close orbital behavior. 

A difference in the behavior of the Hernquist models 
when compared to the King models is the bending of the 
evolutionary tracks at lower masses. The decay of the ap- 
centric radius stops in both models, but the Hernquist model 
dwarfs continue to loose their mass after the orbit evolution 
is terminated. All of the models but HI are completely dis- 
rupted in our simulations during a finite time. 

The central mass distribution of a dwarf can not al- 
ways be represented by the Hernquist profile. The center of 
a dwarf may have a much denser nucleus or even a black 
hole. In this case, the dense central parts of the dwarf may 
survive the disruptive gravitational shocks and settle on a 
low-energy orbit of apocentric radius smaller than the So- 
lar circle for instance, or even settle within the bulge of the 
Galaxy. Our simulations follow the mass evolution of the 
dwarf until the mass of the remnant is about few x 10^ Mq 
. As it is seen from Fig. 1111 the dwarf orbits stop their decay 
when the mass of the remnant is below 10* Mq . Therefore, 
regardless the mass evolution of the dwarf at the very final 
stages, a hypothetical central object will be launched at the 
dwarf orbit when mass of the remnant is about 5 x 10^ Mq 
. The parameters of the final orbit of the central nucleus of 
the dwarf can be estimated rather firmly. In particular, the 
central nucleus of the Hernquist dwarf which has the high 
volume density will have an orbit with apocentric distance 
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Figure 13. Evolution of the angular momentum as a function of 
the decaying apocentric radius of the Hernquist dwarfs. 



of 10 ~ 20 kpc. Remarkably, model H4 sinks to the Galactic 
center with apocentric radius of 6 kpc which is close to an 
apocentric distance of uo Cen. It has to be noticed, however, 
that a progenitor galaxy has to contain in its center a com- 
paract self-gravitating nucleus which can survive disruptive 
encounters with the bulge and the disk of the Milky Way. 

The disk thickening due to an artificial heating, and 
due to an interaction with the satellite can affect the orbit 
decay. This effect depends though on the satellite orbit. In 
our simulations which model the origin of Omega Centauri, 
the orbit decay occurs during first 2 Gyr, or during five-six 
dwarf revolutions around the Milky Way center (Tsuchiya et 
al. 2003). We find that during this time, the satellite orbit 
is not co-planar with the Milky Way disk, and the dwarf 
spends a short time inside the finite thickness disk of the 
Milky Way. In other words, the disk remains 'thin' for the 
dwarf, and the gravitational drag force acting at every point 
of the dwarf orbit does not change much compared to the 
drag force for the disk with the fixed thickness. 

The final apocentric distance of the nucleus is not a 
monotonic function of the initial mass of the dwarf. The 
model with the smallest initial mass HI ends with an apoc- 
entric radius of 18 kpc, while the models with larger initial 
masses end their course closer (H4) as well as further (H8) 
to the Galactic center. 

Fig. ll2l shows the evolution of the bound mass as a func- 
tion of time which is plotted at the apocenters. All models 
except model HI have truncated evolutionary curves with 
the dwarf completely disrupted in a finite computational 
time. 

Figure 13 shows the specific angular momentum as a 
function of apocentric distance of the sinking dwarf. The 
evolution of the specific angular momentum of the Hernquist 
models is not very different from that of the King models. 



4 HERNQUIST MODELS VS KING MODELS 

We see from our simulations a difference in the behavior of 
the King and the Hernquist models. The mass loss rate for 
the King model dwarfs slows down as the bound mass de- 
creases, while for the Hernquist model dwarfs the mass loss 
rate increase. The Hernquist model dwarfs are more fragile 



than the King dwarfs regardless of their mass and size be- 
cause they have shallower central concentration. This, how- 
ever, only partly explains the fragility of the Hernquist mod- 
els. We present in this section simple arguments explaining 
the difference between the King and the Hernquist models. 

There are two factors that affect central density distri- 
bution of the dwarfs. The gravitational shocks inject energy 
into the dwarf's central regions. Another factor is an indirect 
influence of gravitational shocks on the density distribution 
in the dwarf's central regions. The central density distribu- 
tion changes as a response to the loss of the outer layers of 
the dwarf when the dwarf adjusts itself to a new equilibrium. 
Although these two effects act simultaneously, we estimate 
them separately in our simple analytical approach. 

At the later stages of its evolution, the dwarf orbit has 
a low inclination toward the disc and the bulge shocking 
is larger than that from the disc. An approximate expres- 
sion for the energy injection due to a single bulge shocking 
(Spitzer 1958) has the form: 



A AfdwarfrJ^ /2GAfbulgc\ 

^'^= 3 i^^j ^^^^ 



(10) 



Herep, V , rm, and /3 are the pericentric distance (~ 1 kpc), 
the relative velocity (~ 300 km/s), the mean radius of the 
dwarf and the ratio of the duration of the shock to the stellar 
orbital period 



/3 = 



P «rms(r) 



(11) 



V r 

The correction factor to the impulsive approximation, L{P), 
is determined by the relative interaction time /3 such that: 

L(/3) = l, for /3<1; i(/3) < 1, for /3 > 1. (12) 

Equation (10) can be further approximated by replacing the 
dependence on Mdwarf and in eq. 1101 by that of the 
mass of the dwarf within radius r. 

The central density distributions for the King and the 
Hernquist models can be well approximated by a power-law 
distribution 



p{r) = ^0 (^) 



n — —2 or 



This profile is a simplification of our King model (n 



(13) 

-2) 
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and the Hernquist model (n = —1)- Calculating gravita- 
tional potential energy for the density distributions (13), we 
can find a binding energy of the central sphere of radius r 
for the King and the Ifernquist models respectively: 

GM{rf 



10' 



-E'bouiidj^ — 

and 

-E^boundfj — 



2r 



GM{rf 
3r 



(14) 



(15) 



With help of the above expressions, one can write down the 
fractional energy injection for the King model as: 

"■^ GA/buigc 



AE 

bounds 



16 Ml 



bulge 



3 A/t< 



\ r J p 



p3 pY2 



'-L(P), 



(16) 



Here, the mass of the King central sphere of radius r has 
been expressed using total mass Mtot, and the half mass 
radius r^: M{r) — Mtot{r/2rh)- This gives, for example, for 
the King model K4: 



AE 



3.7 



m 



B'boundK \ IkpC ^ 

Similarly, we find for the Hernquist model H4 
AE 



E, 



5.5 



bouiidj^ 



Ikpc 



(17) 



(18) 



For the accepted parameters, the energy injection is smaller 
than the binding energy if r < 0.5 kpc for the King model, 
and if r < 0.2 kpc for the Hernquist model. 

More important effect is a subsequent adjustment of the 
remnant to a new equilibrium. Let's assume that a Hern- 
quist model is instantaneously stripped to a cutoff radius 
Tcut < ro. In this simplified model we consider mass strip- 
ping as the decrease of the cutoff radius rent - Using the Jeans 
equation, we can write the expression determining the veloc- 
ity dispersion of the dwarf at a central region r as a function 
of cutoff radius Vcut'- 



27vGpori (l - rVrLt) 
2nGpororln Vcut/r, 



n = -2, 
n = -1. 



(19) 



For the King profile, the central velocity dispersion does 
not depend much on the cutoff radius r^ut, while for the 
Hernquist profile, the equilibrium velocity dispersion should 
be lower for smaller rcut- The King model dwarfs keep thus 
their central density distribution almost intact. The inner 
regions are close to equilibrium even though the outer lay- 
ers are removed. On the other hand, the Hernquist model 
dwarf expands after the shock because it has an excess of in- 
ternal kinetic energy necessary for equilibrium. This would 
decrease the density of the dwarf, making it more fragile to 
future gravitational shocks which leads to a complete dis- 
ruption of the dwarf in a relatively short time. 



5 DISCUSSION 
5.1 Living halo 

In all of our simulations presented in this paper we assumed 
that the halo has a fixed potential. This assumption does not 
allow us to take into account the drag force due to halo dy- 
namical friction that can change the dynamics of the dwarf. 
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Figure 14. Comparison of the evolutionary histories between 
models with and without the halo responsiveness. In the apoc- 
enter - bound-mass diagram both models take nearly the same 
track. 



For comparison, we ran one simulation with a living halo 
for dwarf model H 4. The halo gravity was calcu lated with 
the SCF algorithm llHernauist and OstrikeJl992l) . The halo 
is sampled by 100,000 particles, and the self-consistent dy- 
namics of the whole s ystem is thus tre ated by a hybrid al- 
gorithm, SCF- TREE llVine and Sieurd sson 1998; Tsuchr^ 
2002). As expected, we find that the orbit decays faster, and 
the mass loss rate is higher when compared to that of the 
models with the fixed halo. For the fixed halo model, the 
apocenter of the H4-dwarf becomes smaller than 10 kpc af- 
ter 3 Gyr, while in the living-halo simulation it happens after 
1.8 Gyr. However, the evolutionary tracks for the fixed and 
living halo models are remarkably similar in the bound mass 
- apocentric distance diagram (Fig. ll4ll . The shock strength 
depends on the pericentric passage velocity, which in turn 
determines the apocentric distance of the dwarf orbit. As a 
result, the bound mass - apocentric distance evolutionary 
tracks are nearly independent of halo model. Since our ba- 
sic conclusions are based on the bound mass - apocentric 
distance diagrams, we correctly predict the final orbit of the 
remnant. 



5.2 uj Centauri 

An important application of our simulations is the under- 
standing of the origin of u Gen. As discussed in Section 4, 
an uj Cen - like object cannot be created from a dwarf with 
a highly concentrated King density profile. On the other 
hand, the dwarfs with Hernquist density profiles exhibit ac- 
celerated mass loss and are completely disrupted in a few 
Gyr. A solution to this dilemma might be the existence of a 
dense compact nucleus at the center of the Hernquist dwarf. 
Such a dense compact nucleus can sink together with the 
disrupting, decaying dwarf and may survive the disruptive 
encounters with the disc and bulge of the Milky Way. A nu- 
cleus with a mass of about lO'^ Mq does not affect the density 
distribution of the dwarf, and most of the evolutionary his- 
tory will be similar to that of a non- nucleated dwarf galaxy . 
This scenario has been examined in lTsuchiva et al] i2003l) . 
We have placed a massive particle of IO^Mq at the center 
of the Hernquist model H4, and followed its self-consistent 



Disruption of a Dwarf Galaxy Under Strong Shocking: The Origin of ou Centauri 11 



evolution. The nucleus settles into an orbit with apoconter 
of about 6 kpc, and has a bound mass approaching IO^Mq 
. A nucleus as a single massive particle is obviously a simpli- 
fication of the model which does not allow us to model real- 
isticallly the dynamics of the stripped nucleus of the dwarf. 
Although, if the N-body nucleus is compact enough, e.g. 
can be represented by a King model with a core radius of 
~ 10 pc, such a nucleus will survive for a few Gyr on its 
present orbit. Nucleated dwarf galaxies are common among 
dwarf galaxies, and such a scenario seems to be a plausible 
explanation of the origin of ui Cen. 

We arbitrarily modeled a compact nucleus by a soft- 
ened particle with a 35 pc half mass radius, equal to the 
gravity softening length. Zhao (2003) noticed recently that 
such a nucleus is more fluffy than the internal central den- 
sity concentration of the Hernquist model, and might be 
destroyed by the tidal forces from the surrounding dwarf. A 
more compact nucleus would obviously remedy the problem. 
However, with our 35 pc resolution, modeling of a more com- 
pact nucleus would not make any difference. Additionally, a 
compact nucleus will not affect an orbital and mass history 
of the satellite. We note also that the half-light radius of ui 
Cen is 4.8' (de Marchi, 1999) or about 7 pc. 

The subsequent evolution of the proto - Omega Cen 
object is determined by a number of processes: dynamical 
friction, evaporation driven by two-body relaxation, gravita^ 
tional shocks, and mass loss driven by stellar evolution (Fall 
& Zhang 2001). Our numerical resolution and the assumed 
model do not allow us to discuss other effects that might 
affect the future evolution of the remaining object. Gnedin 
et al. (2002) estimate the mass loss from u> Cen due to the 
stellar evolution to be 10~^Mci Gyr^^. Gravitational shocks 
which dominate the mass evolution of the massive globular 
clusters also have a time scale of about a few Gyr (Fall & 
Zhang 2001). The dynamical friction is negligible for lj Cen. 
We can speculate therefore that an object which has mass 
close to that of ui Cen and which has a compact enough 
density distribution can survive on its orbit during a few 
Gyr. 

Recently, Bekki and Freeman (2003) have discussed a 
similar model of the origin of u) Cen. They modeled the 
merger of a nucleated dwarf galaxy with an ancient bulgeless 
Galaxy whose thin disk has twenty percent of today's thin- 
disk mass, and included effects of star formation in the dwarf 
dynamics. They find that the Galactic tidal force induces 
radial inflow in the center of the dwarf that triggers rm- 
clear starbursts. Bekki and Freeman (2003) conclude there- 
fore that the effecient dwarf rmclear chemical enrichment 
can be associated with the origin of the observed metal-rich 
stars in u) Cen. 



6 CONCLUSIONS 

We have studied the dynamical evolution of a falling dwarf 
galaxy, which in turn is influenced by the strong gravitar- 
tional shocks from the bulge and disc of the host galaxy. 
We focused our study on the mass-loss history, and on the 
orbital migration of the dwarf. The dwarf galaxy sinks to- 
ward the center of the host galaxy due to dynamical fric- 
tion, which depends on the dwarf's mass. Simultaneously, 
the dwarf is losing mass due to strong bulge and disc gravi- 



tational shocks. As a result, the mass-loss history and orbital 
decay are closely connected and depend in turn on the in- 
ternal density distribution of the dwarf. We have examined 
two different dwarf models, namely a highly-concentrated 
King model and the Hernquist model. The Hernquist model 
has a density cusp while the King model has a density 
profile nearly proportional to r~^. This difference leads to 
qualitatively different evolutionary courses for the two mod- 
els. Regardless of the mass and the size, all the King models 
stop orbital sinking in a few Gyr, once the remnant mass 
has reached ~ 10* Mq . The central regions of the King 
model dwarf remain gravitationally bound, and the rem- 
nant contimies to orbit for a long time without considerable 
mass loss. Therefore, if the initial density distribution of 
the dwarf galaxy is close to isothermal, its remnant can not 
have a mass much smaller than 10* Mq . On the other hand, 
the Hernquist model dwarfs keep losing mass even after dy- 
namical friction becomes ineffective. As a consequence, the 
remnant can be stripped down to a small object. If there 
is a compact gravitationally bound nucleus in the center of 
the dwarf with a mass of < 10^ Mq and will the size of a 
globular cluster, it will survive for few Gyr on the present 
orbit of uj Centauri (e.g.. Fall & Zhang 2001). Our results 
support the capture scenario for the origin of the globular 
cluster a) Centauri. 
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